# points.py
import numpy as np
import matplotlib.pyplot as plt

class Point:
    def __init__(self, position, M0, T2, B0):
        self.position = np.array(position)
        self.M0 = M0
        self.T2 = T2
        self.B0 = B0
        self.amplitude = M0
        self.phase = 0
        self.history = []

    def clear(self):
        self.amplitude = self.M0
        self.phase = 0

    def apply_y_gradient(self, x_gradient ,y_gradient, duration):
        # 计算相位变化量并添加到相位属性中
        phase_shift = y_gradient * self.position[1] * duration + x_gradient * self.position[0] * duration
        base_phase_shift = 0  # 基础磁场引起的相位变化
        self.phase += phase_shift + base_phase_shift
        self.history.append(('x_y_gradient', self.phase))
        decay = np.exp(-duration / self.T2)  # T2 衰减 这里是下降到多少百分比
        self.amplitude = self.amplitude * decay

    def read_out(self, x_gradient, duration, sampling_rate):
        # 计算复信号序列并返回
        t = np.linspace(0, duration, int(sampling_rate * duration))
        decay = np.exp(-t / self.T2)  # T2 衰减
        # decay = 1
        base_phase = 0  # 基础磁场引起的相位变化
        signal = self.amplitude * decay * np.exp(1j * (self.phase + base_phase + x_gradient * self.position[0] * t))
        self.history.append(('read_out', signal))
        return t, signal

    def main(self):
        # 创建一个点对象，假设 T2 = 0.05 秒，B0 = 1.0
        point = Point([1, 1], 1, T2=0.05, B0=1000.0)

        # 测试 clear() 方法
        point.clear()
        print(f"After clear(): amplitude = {point.amplitude}, phase = {point.phase}")

        # 测试 apply_y_gradient() 方法
        point.apply_y_gradient(y_gradient=1.0, duration=0.01,x_gradient=-1.0)
        print(f"After apply_y_gradient(): phase = {point.phase}")

        # 测试 read_out() 方法
        t, signal = point.read_out(x_gradient=1.0, duration=0.5, sampling_rate=1000)
        print(f"Readout signal: {signal[:5]}...")  # 打印前5个信号值以验证

        # 画出信号的实部和虚部
        plt.figure(figsize=(12, 6))

        plt.subplot(3, 1, 1)
        plt.plot(t, signal.real)
        plt.title('Real Part of Signal')
        plt.xlabel('Time (s)')
        plt.ylabel('Amplitude')

        plt.subplot(3, 1, 2)
        plt.plot(t, signal.imag)
        plt.title('Imaginary Part of Signal')
        plt.xlabel('Time (s)')
        plt.ylabel('Amplitude')

        plt.subplot(3, 1, 3)
        plt.plot(t, np.abs(signal))
        plt.title('Magnitude of Signal')
        plt.xlabel('Time (s)')
        plt.ylabel('Amplitude')

        plt.tight_layout()
        plt.show()

    if __name__ == "__main__":
        main()